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SUMMARY 

m,.ltmi^n^ Pr n Pe l la ? t CFD dC u ign COde haS been devel °P ed to Predict the local mixing of 
multiple propellant streams as they are injected into a rocket motor. The tripropellant system 

of primary interest is G0 2 /H 2 /RP-1. The code utilizes real fluid properties to account for the 
mixing and finite-rate combustion processes which occur near an injector faceplate. Thus the 
analysis serves as a multi-phase homogeneous spray combustion model, which explicitly avoids 
the droplet tracking analysis that makes alternative analyses extremely slow to compute. Of 

pm U nh!^ SimP1 I 1Ca !i° n < J 0Uld hJ ! Ve been affected hy compromising the combustion description and 
p asizmg the droplet tracking. However, proper accounting of the combustion allows 
accurate gas-side temperature predictions which are essential for accurate wall heating analyses 
The complex secondary flows which are predicted to occur near a faceplate cannot be 
quantitatively predicted by less accurate methodology. Since very few tripropellant injector 
xpenments have been reported, the methodology has also been applied to bipropellant injectors 

‘ vahd . atl0n - Test cases have been simulated to describe an axisymmetric 
tnpropellant coaxial injector and a 3-dimensional RP-l/LO, impinger injector system. The 

analysis has been shown to realistically describe such injector combustion flowfields, but further 
exerasing of the code is necessary to accomplish its computational optimization. Since the 
® tbe .P red ^ ted combustor flowfields exceed the state-of-the-art capabilities to measure 
. flowfields, the code is also valuable to design meaningful future experiments by 
determining the critical location and type of measurements needed. 
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DEVELOPMENT OF TRIPROPELLANT CFD DESIGN CODE 


1.0 INTRODUCTION 

SECA, Inc. has developed a CFD design code to describe the operation of a tripropellant 
injector system, which, for example, utilizes both RP-1 and H 2 as fuels. This model is 
sufficiently general to define the mixing characteristics of several injector concepts. Future 
engine designs may be devised with the code developed in this investigation. The code provides 
the local mixture ratio and state of the propellants near the injector face so that the remainder 
of the combustor flow can be calculated with conventional multi-species gas-only CFD codes. 
It is anticipated that laboratory test data will become available to validate the CFD injector 
model, although adequate data do not currently exist. 

The objective of developing an advanced tripropellant CFD design code is to provide an 
analytical method to predict the local mixing of the multiple propellant streams as they are 
injected into the motor. This analysis must be described with a computer code which is 
sufficiently fast to provide very detailed predictions and structured to be used as a submodel in 
a general three-dimensional CFD simulation of the entire motor operation. It is recognized that 
other spray combustion models have been and are being developed, but such models have not 
matured to the point that efficient rocket performance and heating analyses can be made for 
practical rocket motor designs. Extensive data bases and much experience are available for 
designing bipropellant rocket motor injectors, but such valuable information is not available for 
tripropellant and advanced bipropellant injectors. Furthermore, the new CFD design code can 
be used to reduce the large number of costly experimental injector studies which would 
otherwise be required to design future rocket motors. 

To describe tripropellant combustors, real fluid properties representing both sub- and 
super-critical conditions for RP-1, H 2 , and 0 2 must be used. A combustion kinetics model for 
RP-1 including soot formation and oxidation was developed. Radiation from sooty combustion 
gases is important and was assessed, but, due to a reduction in program funding, was not 
included in the final CFD model. To treat fluids in arbitrary states, the pressure correction step 
in the FDNS CFD code was completely changed from the method which has been used for ideal 
gases and incompressible liquids. The tripropellant CFD design code resulting from this 
development is termed the Finite Difference Navier-Stokes - Real Fluids Version (FDNS-RFV) 
code. 


The tripropellant design code produced from this research meets the objectives of the 
proposed study and may be termed a homogeneous spray combustion code. Since the mixed 
gas/liquid flowfield is analyzed without identifying and tracking drops, the gas and liquid are 
always assumed to be in thermal equilibrium. This assumption greatly simplifies the calculation 
while still approximating the multi-phase flow near the injector face. The limitations and 
alternatives to using this assumption are discussed in this report. Both shear coaxial and 
impinging doublet injectors were simulated with the code as demonstration cases. 
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2.0 METHODOLOGY 

During the course of the Space Shuttle Main Engine (SSME) and the Fastrack engine 
design processes, severe erosion problems on the liner wall near the injector face have occurred. 
Study and understanding of the spray combustion flows have posed a great challenge to the 
injector, and hence, the rocket engine design. The effect of the spray vaporization is also 
important in analyzing the internal ballistics of rocket motors. Though extensive experimental 
studies can help designers/researchers understand and analyze the flow characteristics, building 
various injector/motor hardware can be very costly. Hence, an efficient CFD spray combustion 
code can improve the rocket engine design and reduce the cost of hardware. The objective of 
the present work is to develop an efficient CFD design tool to simulate real fluid injector 
flowfields for the liquid and hybrid rocket engines. Lately, CFD methodology has improved to 
the point that numerical analyses of turbulent, reacting flows provide practical design models for 
combustion devices. However, multi-phase (liquid-spray) combustion CFD models are still 
under development. The main reason is that not many experimental studies for spray 
combustion flows have been conducted, and, unfortunately, most of those tests did not provide 
detailed flow data. Hence, the underlying physics is still not well understood, and this 
shortcoming is reflected in the available CFD spray combustion models. Due to the importance 
of the benchmark test data for spray combustion flows, NASA Marshall Space Flight Center 
funded Pennsylvania State University to conduct experimental tests to provide flowfield 
measurements of the gas/liquid coaxial injector 1 ' 2 . The availability of such data allow the present 
homogeneous CFD spray combustion model to be benchmarked to the extent that sufficient 
measurements were made. Deleting the droplet tracking simplifies the CFD code to the point 
where meaningful 3-D flowfield simulations and design studies can be made. 

Tripropellant injection may be accomplished with elements comprised of like- or unlike- 
doublets, triplets, or other variations of coaxial injectors. All of these configurations require that 
three-dimensional CFD analyses be performed to describe the initial mixing of a single injector 
element. The entire injector assembly must then be described with an ensemble of such injectors 
or with alternative methodology. 

The complete combustor/prebumer tripropellant CFD design code (FDNS-RFV) will 
require: (1) definition of the inlet boundary conditions to each representative injector element, 
as obtained from a model of the feed system such as that provided by SECA’s porosity model 
for the SSME main combustion chamber 3 , (2) the CFD analysis of each representative injector 
element (the number of representative elements should be small, for a given motor design); and 
finally (3) the use of the injector model results as boundary conditions for a CFD simulation of 
the entire motor and nozzle flow. The injector element analysis shall be accomplished with a 
CFD calculation. The injector element analyses shall provide streamtube bundles to initiate a 
CFD simulation of the entire combustor. Both the injector element analyses and the entire motor 
analysis can be used to determine local wall heating effects. 

2.1 Real Fluids Version of FDNS (FDNS-RFV) 

The framework of the CFD model is an elliptic, Finite Difference Navier-Stokes (FDNS) 
flow solver 4 ' 5 , which employs a predictor plus multi-corrector pressure-based solution algorithm 
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so that both compressible and incompressible flow problems can be analyzed. High order 
upwind, total variation diminishing (TVD), or central difference schemes plus adaptive second- 
order and fourth-order dissipation terms are used to approximate the convection terms of the 
transport equations such that computed supersonic flows through a shock regime will be 
stabilized. Second-order central differencing schemes are used for the viscous and source terms 
of the governing equations. Vectorized point implicit, conjugate gradient, and generalized 
minimal residual matrix solvers 6 (GMRES) are optionally employed to insure a stable, accurate, 
and fast convergence rate. Multi-block, multi-zone options are included in the FDNS code so 
that problems with complex geometries can be analyzed efficiently. An extended k-e turbulence 
model with a modified near-wall boundary treatment 7 is utilized in the code for turbulent flow 
computations. In this approach, the profile of the non-dimensionalized velocity u + is formulated 
based on the velocity profile suggested by Liakopoulos 8 . This velocity profile provides a smooth 
transition between logarithmic law of the wall and linear viscous sublayer variations. In 
addition, a modified wall function treatment 9 , which is suitable for compressible flows by 
accounting for the aerodynamic heating effects, had previously been incorporated in FDNS to 
provide appropriate heat transfer calculations near wall boundaries. To simulate spray 
combustion flows, several additional submodels, as described in subsequent paragraphs, were 
developed and incorporated into the FDNS code under the Phase II SBIR investigation 10 which 
this Phase III study is derived from. The Phase II 10 progress reports also provide more detail 
for the methodology reported herein. 


In order to analyze the effect of phase changes in the spray combustion flows, real-fluid 
thermal and caloric equations of state were developed 10 . Interfacial heat and mass transfer 
effects were neglected for the very high pressure (near or above critical pressure) in existing 
rocket engine combustion chambers. The HBMS equations of state 1112 were selected for this 
purpose. These equations are: 
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where P c , T c , p c , and Z c are the pressure, temperature, density, and compressibility at the critical 
condition. By are the coefficients of the thermal property polynomial for a given species. H and 
Ho are the real-fluid enthalpy and ideal-gas enthalpy for a given species, and R is the gas 
constant. These equations are not only of sufficiently high order that properties are accurately 
predicted for a wide range of conditions, but component submodels may be easily modified. In 
this instance, the vapor pressure curve and the liquid phase density correlations were improved 
over the original HBMS formulation. Multicomponent mixtures were treated by adding partial 
specific volumes or pressure. The partial volume methodology is essential to provide accurate 
and robust predictions when a small amount of multi-component vapor and a large amount of 
liquid, or vice versa, are present at the same location. The properties routines also include 
correlations for the transport properties. Furthermore, this real-fluid property submodel can be 
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ported into other CFD code, or used as a stand alone program. 


Most of the pressure-based methods are constructed based on the perfect gas law, which 
is no longer valid in simulating the spray combustion flows. In the present model, flow mixture 
density is calculated in two steps. In the predictor step, a constant pressure procedure is 
employed to solve flow velocity components, mixture energy and species concentrations, and 
then the flow mixture density and temperature are calculated from the real-fluid model. Like 
other pressure-based methodology, a correlation between pressure and density is used to convert 
the continuity equation to the pressure-correction equation as the corrector step. However, the 
perfect gas law can not be used as the correlation. In the present model, a constant-temperature 
derivative is used to correlate the pressure change and the density: 


P' 

dp ' 

= P'7 

dp 


dP 

T 

dP 


0pP' 


(3) 


where for the perfect gas: 



(4) 


and for the real -fluid model: jS p = (5) 

a 2 

where 7, R and a are the specific heat, gas constant and sonic speed of the flow mixture, 
respectively. Hence, the pressure-correction equation can be derived from the continuity 
equation based on the above correlation, and can be written as 

• ( V/3 p P' ) - V • (p* D p V P ; ) = -V • (p* V) - p 
At At 

where D P is the matrix coefficient of the momentum equations, and the superscripts of * and n 
denote the values at the intermediate and previous time steps. 


2.2 Fluid Properties 

Methods for estimating real fluid properties of L0 2 /RP-1/H2 mixtures and their 
combustion products have been developed 10,13 and included in the FDNS-RFV code. The 
properties of a RP-1 surrogate fuel are presented in Table 1. The Core Labs report on RP-1 10 
approximated the constituents to be about 42% liquid volume of paraffins, 25% 
monocycloparaffins, 24% di- and tri-cycloparaffins, 3% benzenes and 6% naphthalenes. 
Examination of the constituents revealed about 12% liquid volume of molecules with 10 carbon 
atoms (CIO’s), 13% Cll’s, 29% C12’s, 25% C13’s, 11% C14’s and the rest either larger or 
smaller. The surrogate fuel was designed to incorporate these breakdowns into the formulation. 

In order to simplify the surrogate model, the carbon breakdown was reduced to 20% 
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Table 1. Surrogates for RP-1 


Species . 

Type 

Formula 

Vol % 

Den 

Mole Fr 

n-UNDECANE 

Paraffin 

C11H24 

4.70 

.781 


DODECANE 

Paraffin 

C12H26 

6.00 

.791 


n-TRIDECANE 

Paraffin 

C13H28 

18.80 

.819 

.17828 

n-TETRADECANE 

Paraffin 

C14H30 

12.50 

.761 


n-HEXYLCYCLOPENTANE 

Monocyclic Paraffin 

C11H22 

2.70 

.782 

.02921 

n-HEPTYLCYCLOPENTANE 

Monocyclic Paraffin 

C12H24 

3.60 

.782 

.03570 

n-OCTYLCYCLOPENTANE 

Monocyclic Paraffin 

C13H26 

11.20 

.796 

.10437 

n-NONYLCYCLOPENTANE 

Monocyclic Paraffin 

C14H28 

7.50 

.803 

.06547 

BICYCLOPARAFFIN 1 

Polycyclic Paraffin 

C11H20 

11.30 

.852 

.13496 

BICYCLOPARAFFIN2 

Polycyclic Paraffin 

C12H22 

14.70 

.819 

.15453 

PENTAMETHYLBENZENE 

Mononuclear 

C11H16 

1.30 

.806 

.01509 


Aromatic 





HEXAMETHYLBENZENE 

Mononuclear 

C12H18 

1.70 

.786 

.01758 


Aromatic 





DIMETH YLN APHTH ALENE 

Dinuclear Aromatic 

C12H12 

4.00 

.967 

.05285 

Formula 



C12.39H24.15 


H/C 



1.9492 



Molecular Weight 



173.1453 


| Liquid Density 



.8111 gm/cc 


Breakdown of Number of Carbon Atoms per Molecule 





Fraction of LV of Cl Is (.20) 



20.0000 



Fraction of LV of Cl 2s (.30) 



30.0000 



Fraction of LV of Cl 3s (.30) 



30.0000 



Fraction of LV of C14s (.20) 



20.0000 



Breakdown of Liquid Volume by Type 





Paraffins (.42) 



42.0000 



Monocyclic Paraffins (.25) 



25.0000 



Polycyclic Paraffins (.24) 



26.0000 



Mononuclear aromatics (.03) 



3.0000 



Dinuclear aromatics (.06) 



4.0000 



Stoichiometric Products: 






Moles of 0 2 burned 



18.426 



Moles of C0 2 formed 



12.389 



Moles of H 2 0 formed 



12.075 




Heat of Formation of Surrogates -.56940E+02 kcal/mol 

-.32886E+02 kcal/lOOgm 

Heat of Formation of Products 1863 IE +04 kcal/mol 

-.10760E+04 kcal/lOOgm 

Heat of Combustion 1 806 1 E + 04 kcal/mol 

-.10431E+04 kcal/lOOgm 
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Cll’s, 30% C12’s, 30% C13’s and 20% C14’s. Eighteen candidate constituents were identified 
which satisfied both Core Lab breakdowns. Various liquid volume fractions of each were 
analyzed and compared to estimated RP-1 molecular weight (175), hydrogen to carbon ratio 
(H/C = 1.95) and liquid density (0.83 grams/cc). A reasonable value of H/C was unattainable 
due to the low H/C’s for aromatics. Since the reported percentage (9%) of aromatics was higher 
than reported by other sources, and since complex dinuclear aromatics (naphthalenes) are often 
indistinguishable from polycyclic paraffins, the naphthalenes percentage was reduced by 2% and 
the polycyclic paraffin percentage was increased by 2% . This resulted in an H/C close to 1.95. 
These considerations are reflected in the surrogate shown in Table 2. 


Table 2. Critical Properties for RP-1 Surrogate 

Species 

Mol Fr. 

Pc(bar) 

Tc(°K) 

Vc(gmol/cc) 

Slilllill 

n-UNDECANE 

.05013 

19.70 

638.80 

647.05 

.240 

DODECANE 

.05948 

18.20 

658.20 

694.58 

.231 

n-TRIDECANE 

.17828 

17.20 

676.00 

784.25 

.240 

n-TETRADECANE 

.10235 

14.40 

693.00 

920.29 

.230 

n-HEXYLCYCLOPENTANE 

.02921 

21.30 

660.10 

590.05 

.229 

n-HEPTYLCYCLOPENTANE 

.03570 

19.40 

679.00 

689.67 

.237 

n-OCTYLCYCLOPENTANE 

.10437 

17.90 

694.00 

744.64 

.231 

n-NONYLCYCLOPENTANE 

.06547 

16.50 

710.50 

805.54 

.225 

BICYCLOPARAFFIN 1 

.13496 

27.79 

712.20 

539.09 

.253 

BICYCLOP ARAFFIN2 

.15453 

24.36 

721.80 

593.72 

.241 

PENTAMETHYLBENZENE 

.01509 

24.85 

719.00 

541.27 

.225 

HEXAMETHYLBENZENE 

.01758 

22.40 

758.00 

607.72 

.216 

DIMETHYLNAPHTHALENE 

.05285 

31.09 

771.00 

517.53 

.251 

Mixture Critical Temperature (K) 



692.776 



Mixture Critical Pressure (bar) 



19.974 



Mixture Critical Volume (cc/gmol) 


646.093 



Mixture Critical Zc 



.224 




The true critical properties for the RP-1 surrogate fuel were calculated using Gibbs 
method 13 . Unknown saturated liquid volumes of the constituents were estimated using the 
method of Gunn and Yamada 14 . Unknown critical properties were estimated using the group 
contribution method of Joback 15 . As before, Gibbs method resulted in a high value for the 
mixture critical volume. The literature 1416 suggest using the resultant critical pressure and 
temperature but not the critical volume when Gibbs method is employed with cubic equations 
of state. An alternative was to use the method of Chueh and Prausnitz 14 to estimate critical 
volume. This method utilizes an interaction parameter which varies depending on whether the 
constituent is a paraffin, cycloparaffin, aromatic, etc. These interaction parameters cause the 
critical volume to vary wildly, hence care must be used in evaluating them. For the complex 
molecules, such as n-nonylcyclopentane or hexamethylbenzene, more than one category may be 
applicable. These methods were utilized and the interaction parameters adjusted until a 
reasonable value for the critical volume was obtained. The results are presented in Table 2. 
A plot of the resulting thermal equation of state is shown in Figure 1. 
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Figure 1 

Pressure and Density Along Isotherms 

RP — 1 Comprised of 1 3 Surrogates 


2.0 1.4 1.20 1.15 
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2.3 Combustion Kinetics 

Combustion kinetics rates for H 2 /RP-1 at high pressures are very fast except possibly for 
the formation and oxidation of soot and for the pyrolysis of RP-1 constituents. If local chemical 
equilibrium were assumed, ignition delays upon injection could not be predicted, no soot would 
be formed if any oxygen were present, and all large molecular species in the RP-1 would be 
immediately pyrolyzed. The RP-1 contains a mixture of so many species that they cannot be 
rigorously identified. Kinetics models for RP-1 oxidation that treat the hydrocarbons as either 
large aromatics, large aliphatics, or small two-carbon fragments are desirable, but such models 
are too complex for direct use in a CFD code. The species which comprise RP-1 change from 
lot-to-lot of the fuel, but these composition changes are unknown, therefore the one surrogate 
fuel characterization was used. One global pyrolysis/partial-oxidation step for the surrogate fuel 
was used as the RP-1 kinetics model. The rate used for this global step is given in Table 3. 
The combustion is completed by a wet CO equilibrium model. This RP-1 kinetics model 
accounts for the large amounts of H 2 which would be present in a tripropellant with the 
equilibrium H 2 /0 2 steps in the wet CO reactions. H 2 /0 2 combustion for gaseous 0 2 is modeled 
as an equilibrium process. For the LO z /H 2 cases analyzed in this report, the global reaction 
shown in Table 4 was used. The approximation introduced with this global model is that the 
distribution of combustion products is evaluated for equilibrium combustion at one mixture ratio. 
If the stoichiometric coefficients were expressed as a function of local mixture ratio, the global 
equation simulation could be made quite accurately. However, as it was used, this model is an 
approximation which is considerably better than assuming complete combustion and is 
computationally efficient. 

2.4 Radiation From Sooty Flames 

RP-1 and hybrid fuels tend to form soot during the combustion process. Soot, as well 
as the gases C0 2 , H 2 0, and CO, produce significant thermal radiation. Since the soot particles 
are expected to be small, the soot radiation is continuum and is exchanged by emission and 
absorption only, as is the radiation for the gases. However, gaseous radiation is strongly and 
discontinuously wavelength dependent. To predict radiation in combustion devices, three factors 
must be considered: (1) radiation properties for the optically active species, (2) the method of 
calculating the radiative exchange, and (3) the accuracy of the concentration predictions of the 
radiating gas species and particles. The essential elements which define these factors are 
presented in Appendix A. No attempt has been made to compromise the radiation analysis 
presented in Appendix A solely for computational expediency. 

The critical issue is the accurate prediction of soot and liquid masses and drops. 
Absolutely no data to indicate the magnitude of the radiation shielding effects of liquid particles 
in the vicinity of the injector face-plate are available. Gaseous species must also be predicted, 
but these predictions are both more accurate and less critical since gases account for less of the 
radiation. The flowfield analyses described in this report are useful for making such predictions; 
however, identification of typical combustion systems and validation data must be accomplished. 
The CFD solutions are too cumbersome to yield the necessary flowfield data for an 
indiscriminate parametric investigation of combustors. Definitive experiments to define radiation 
in tripropellant rocket combustion devices are non-existent. Since soot predictions and their 
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verification are the final result of tripropellant flowfield analyses and since this investigation was 
foreshortened because of funding reductions, the radiation methodology described in Appendix 
A was not used to make radiation predictions for tripropellant systems. The radiation aspects 
of tripropellant systems should be revisited when (and if) more firm design concepts are 
established. 


Table 3. Combustion Kinetics Models for RP-l/O, 


RP-1 + 6.195 0 2 - 12.39 CO + 12.075 H 2 
Rate= AT B exp{-E/RT} [RP-1] 05 [OJ 10 

A = 4.4963xl0 9 , B = 1, E/R = 2.679X10 4 


RP-1 - 12.39 C. + 12.075 H 2 

Rate= A e*/>{-25(l-T/T 0 ) 2 } [RP-1] 143 [0 2 + e] 05 

A = 0.5, T 0 = 1700°K, e = O.OOlp/M^ g-mole/cm 3 


C. + 0.5 0 2 -> CO 

6 M„ RT f K ,rlf 

Rate= w> _ - -j— — p— + K 2 (1 - \p) [C.][0 2 ] 

p, U, A + *4 F 0j 


P o K 3 

where \p = 1 + — - — 

*2 


; and K ; = A^ 1 ^ , i = 1 -* 4 


A, = 20, E,/R = 1.509X10 4 

A 2 = 4.46xl0' 3 , E 2 /R = 7.6497x10 s 

A 3 = 1.51x10 s , Ej/R = 4.8817X10 4 

A 4 = 21.3, E 4 /R = -2. 0634x10 s 

R u = 82.06 atm • cm 3 /g-mole • °K 

p, (density of soot) = 1.86 g/cm 3 

D, (particle diameter of soot) = 40 nm (assumed value) 

M Wi (molecular weight of soot) = 12 g/g-mole 

Poa (partial pressure of 0 2 , in atm) = P • (N^: mole fraction of 0 2 ) 


Table 4. Combustion Global Kinetic Model for H 2 /0 2 

Reaction 

A 

B 

E/R 

1.827 H 2 + 0 2 - 1.552 H 2 0 + 0.094 O + 0.196 H 
+ 0.354 OH 

3.1000E22 

0 

2.4070E4 
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3.0 APPLICATIONS TO INJECTOR DESIGN 

At the inception of this investigation, strong consideration was being given to the use of 
tripropellant rocket motors which utilized H 2 , RP-1, and 0 2 . Both impinging and coaxial 
elements were being evaluated. Subsequent NASA programs dropped the tripropellant concept. 
The only tripropellant injector design which has been experimentally investigated was the 
effervescent H2/RP-1 coaxial injector 17 . Therefore, due to the lack of validation data, several 
bipropellant injectors of the type used in the Fastrack motor and bipropellant coaxial elements, 
as well as, the aforementioned coaxial tripropellant injector were simulated to evaluate the 
tripropellant CFD design code. The results of these simulations are reported in this section. 

3.1 The Entire Fastrack Motor 

To initiate the study of clustering injector elements to represent an entire combustor, a 
series of a 2-D axisymmetric numerical global simulations of the Fastrack test motor were 
conducted. The major goal of this study was to determine the temperature along the motor wall, 
in order that wall heating could be subsequently evaluated. The motor configuration is shown 
in Figure 2. The flow conditions which were used to parametrically represent the injector flow 
are shown in Tables 5-7. Uniform flow of combustion gas through most of the injector and 
center flow and film cooling flow of RP-1 were used to represent the injector configuration. At 
first, an ideal gas model was employed to simulate RP-1 and combustion gas. A large gap for 
the center fuel flow and outer film cooling flow was used as shown in Table 5, to compensate 
for the low density RP-1 vapor. For the ideal gas model, two numerical calculations were 
conducted, one was with frozen chemistry and the other with finite-rate chemistry. The results 
shown in Fig. 3, indicate that the difference is very small because the major part of the flow is 
in the chemical equilibrium and is uneffected by the small amount of fuel which flows through 
the center and outer coolant streams. Therefore, the frozen chemistry at equilibrium chamber 
conditions was chosen for further study. Figure 4 shows the results of the ideal gas model with 
frozen chemistry simulation. 

In order to study the effect of a liquid phase for the center fuel and the outer coolant fuel 
flow, the real fluids model was employed to simulate these two streams, as opposed to the ideal 
gas in the previous calculation. Intially, the film coolant was assumed to occupy the space 
between a point midway radially between the film coolant holes and the outer row of primary 
fuel holes. The small-gap case assumed that the film coolant was injected through an area 
defined by a ring of thickness equal to the gap between the film coolant holes and the wall. As 
a result, the predicted wall temperature is much lower than in the previous case and is thought 
to be a better simulation of actual conditions. 

Two different mixing models, the extended k-e model with and without a temperature 
correction were used to study the effect of compressibility on the species mixing and thus on the 
nozzle wall temperature. The temperature correction model modifies the source term S £ , in the 
£ equation as 
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Figure 2 

SKETCH OF FAST TRACK NOZZLE GEOMETRY 
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Table 5. INLET CONDITIONS FOR FAST TRACK NOZZLE SIMULATION 

(IDEAL GAS, LARGE GAP) 


RP-1 flow rate 


0 2 flow rate 


Starting/Ending 
radial location 


Temperature 


Density 


Axial Velocity 


Inner Fuel Region Hot Gas Mixture Outer Fuel Region 


0.76678 lb/sec 


0/0.65" 


532 °R 


15.457 lb/ft 3 


5.382 ft/ sec 


21.16322 lb/sec 


55.78 lb/sec 


0.6574.323" 


6714 °R 


0.16457 lb/ft 3 


1173.25 ft/sec 


3.32 lb/sec 


4.32374.6" 


532 °R 


15.457 lb/ft 3 


3.983 ft/sec 


Table 6. INLET CONDITIONS FOR FAST TRACK NOZZLE SIMULATION 
(REAL FLUID MODEL, LARGE GAP) 


RP-1 flow rate 


0 2 flow rate 


Starting/Ending 
radial location 


Temperature 


Density 


Axial Velocity 


Inner Fuel Region Hot Gas Mixture Outer Fuel Region 


0.76678 lb/sec 


0/0.65" 


532 °R 


54.04 lb/ft 3 


1.539 ft/sec 


21.16322 lb/sec 


55.78 lb/sec 


0.6574.323" 


6714 °R 


0.16457 lb/ft 3 


1173.25 ft/sec 


3.32 lb/sec 


4.32374.6" 


532 °R 


54.04 lb/ft 3 


1.139 ft/ sec 
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Table 7. INLET CONDITIONS FOR FAST TRACK NOZZLE SIMULATION 
(REAL FLUID MODEL, SMALL GAP) 


Inner Fuel Region Hot Gas Mixture Outer Fuel Region 


RP-1 flow rate 


0 2 flow rate 


Starting/Ending radial 
location 


Temperature 


Density 


Axial Velocity 


0.76678 lb/sec 


0 


0/0.45" 


532 °R 


54.04 lb/ft 3 


3.212 ft/sec 


21.16322 lb/sec 


55.78 lb/sec 


0.4574.5" 


6714 °R 


0.16457 lb/ft 3 


1068.99 ft/sec 


3.32 lb/sec 


0 


4.574.6" 


532 °R 


54.04 lb/ft 3 


3.095 ft/sec 


Table 8. NUMERICAL RESULTS OF THE NOZZLE EFFICIENCY 
AND WALL TEMPERATURE FOR VARIOUS MODELS 



CEC Calculation 


Gas RP-1 
(large gap) 


Liquid RP-1 
(large gap) 


Liquid RP-1 
(small gap) 


frozen chemistry 

5474 

5663 

0.946 

3618 

finite rate 

5505 

5705 

0.953 

3750 

extended k-e 

5300 

5483 

0.916 

4158 

T-corrected k-e 

5275 

5457 

0.911 

3742 

extended k-e 

5340 

5555 

0.928 

4334 

T-corrected k-e 

5347 

5563 

0.929 

3550 
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S - P [ca| - CtJ ♦ c, (T-f- <*> 

where C,, C 2 , C 3 , and C 4 are modeling constants, and T* is a ratio of local to ambient 
temperature which when raised to the C 4 power 4 . This model reduces the apparent turbulent 
viscosity across temperature gradients. The results, as shown in Figures 5 and 6, demonstrate 
that with the temperature correction model the film cooling is predicted to be more effective. 
However, compared to the result of the ideal gas model, the liquid RP-1 coolant is less effective 
than the gas RP-1 coolant. This is because the momentum of the liquid RP-1 is smaller than that 
of gas RP-1. In order to model the liquid RP-1 coolant flow more accurately, a small-gap case 
was simulated where the RP-1 coolant is injected at the coolant hole location as described in 
Table 7. As with the previous case, the same two mixing models were employed, and the 
results are shown in Figures 7 and 8. 

Figure 9 and Table 8 summarize these simulations. The lowest predicted wall 
temperatures are believed to be the most accurate. The C* efficiency values are thought to be 
realistic. Additional work is required to evaluate the effect of the ablative liner response to this 
predicted wall temperature distribution. 

3.2 G0 2 /GH 2 Shear Coaxial Injector 

Due to the lack of definitive test data on a tripropellant injector element, several injector 
configurations were simulated to evaluate the homogeneous spray combustion model in a piece- 
meal fashion. Santoro and his co-workers at PSU investigated several coaxial injector 
elements 1,2,17 . These consisted of G0 2 /GH 2 , L0 2 /GH 2 , and RP-1/GH 2 /G0 2 injectors. All of 
these injectors were numerically simulated to evaluate various aspects of the CFD modeling. 

The G0 2 /GH 2 shear coaxial injector test case 1 , was simulated to validate the capability 
of the proposed model to predict mixing and chemical reaction phenomena. This test case is the 
first experiment which was designed to provide very detailed injector flow data. Measurements 
included: mean flow, root-mean-square fluctuation velocities, OH radiance, and stable species 
concentrations for high pressure combustion. A uni-element shear coaxial injector, as sketched 
in Figure 10, was used to introduce the propellants into the combustion chamber. The numerical 
simulations are for a computational domain which starts 25.4 mm upstream and stops 8 inches 
downstream from the injector face. The inlet conditions of the propellants are listed in Table 
9. The H 2 /0 2 combustion was assumed to occur under conditions of local chemical equilibrium 
which are appropriate at these pressures. The numerical results are compared to the test case 
as shown in Figures 11-12. The simulated flowfield is plotted in terms of temperature and is 
illustrated in Figure 13. The agreement between the numerical results and the test data are very 
good except for the H 2 0 mole fraction. In Fig. 11, the growth of shear layer is shown to be 
well predicted which indicates the turbulent mixing and the chemical reaction are simulated 
accurately. The major part of the discrepancy between the measurements and the predictions 
for the G0 2 /GH 2 injector are attributed to inaccuracies in the measurements. For example, if 
the reported water measurements are correct, very little combustion would have occurred. This 
is not consistent with the good agreement obtained between the predicted and measured velocity 
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Figure 10. The Schematic of the Uni-element Shear Coaxial Injector Test Facility 
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Table 9. The Inlet Flow Conditions of the G0 2 /GH 2 Injector Test Case 



go 2 

gh 2 


Pressure (MPa) 

1.31 

1.31 


Temperature (°K) 

290 

298 


Mass flow rate (kg/s) 

0.04196 

0.0104 



Axial velocity (m/s) 

51.1 

177.2 



_ 

Table 10. The Inlet Flow Conditions of the L0 2 /GH 

2 Injector Test Case #1 



lo 2 

gh 2 

— 

Pressure (psi) 

453 

453 


Temperature (°R) 

210 

557 

— 

Density (lb/ft 3 ) 

62.79 

0.1516 


Mass flow rate (lb/s) 

0.397 

0.049 


Table 11. The Inlet Flow Conditions of the L0 2 /GH 2 Injector Test Case #2 


lo 2 

gh 2 

Pressure (psi) 

260 

260 

Temperature (°R) 

210 

557 

Density (lb/ft 3 ) 

62.42 

0.0864 

Mass flow rate (lb/s) 

0.371 

0.073 


profiles. Furthermore, the summation of the measured mole fractions is not close to unity. This 
difference cannot be solely attributed to unmeasured radicals. Temperature measurements were 
not reported, but verbal reports indicated that temperature was estimated to be about 3000 °K 
in the shear layer. The discrepancy between the measurements and the numerical predictions 
for the velocity near the chamber wall is caused by the nitrogen purge used to protect the optical 
windows. The purge flow was not accounted for in the numerical simulations. As far as the 
fluctuating velocity is concerned, the disagreement in the comparison is anticipated because the 
values of the axial fluctuating velocity are calculated from the turbulence kinetic energy, 
assuming isotropy. This assumption is widely recognized to be inaccurate for turbulent jet 
flows, hence its use probably accounts for the predicted isotropic axial fluctuating velocities 
being smaller than the measured values. None of the comparisons made in this evaluation 
suggest that any modifications to the CFD model are necessary. 
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(a) x = 25.4 mm 




(b) x = 50.8 mm 



(c) x = 127 mm 

Figure 11 Comparison Between Measured and Calculated Mean and Root-Mean-Square 
Axial Velocities for Three Axial Locations (equilibrium chemistry) 
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It Moh Fraction 


foU* Fraction 

(a) x = 25.4 mm 


to Wole fraction 
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to 


-10 



ft Fraction 


03 Mole Fraction 

(b) x = 50.8 mm 


to Mob Fraction 





02 Mote Fraction to bob fraction 

(c) x = 127 mm 

Figure 12 Comparison Between Measured and Calculated Species Concentrations for Three 
Axial Locations (equilibrium chemistry) 
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Figure 13 Predicted Temperature Contours for the G0 2 /GH 2 Inject Flow 
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Figure 14 Predicted Temperature Contours for the LOX/GH 2 Injector Test Case #\ 
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Figure 15 Predicted Temperature Contours for the LOX/GH 2 Injector Test Case # 2 
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3.3 L0 2 /GH 2 Shear Coaxial Injector 

Two LO 2 /GH 2 shear coaxial injector flow test cases, also examined by Santoro and his 
coworkers 2 , were predicted to further validate the real fluid homogeneous spray combustion 
model. The expectation was that the homogeneous spray combustion model would underpredict 
the mixing, since droplet generation and dispersion was not simulated. The configuration of the 
experiment was similar to that of the gas/gas injector test case, except that liquid oxygen flowed 
through the inner tube. The inner tube diameter was 3.42 mm, and the gas hydrogen flowed 
through the outer annulus. The inner and outer diameters of the annulus were 4. 19 mm and 
7.11 mm, respectively. The nozzle was changed to provide the desired test pressure. The inlet 
conditions of the propellants for both test cases are listed in Tables 10-11. A global kinetics 
chemistry for H 2 /0 2 , listed in Table 4, was employed to model the spray combustion because 
the cold temperatures suggested that finite-rate effects could be important. The calculated results 
for both cases are plotted as shown in Figures 14-15. 

Test data for sub-critical combustion were reported to indicate a growing jet of LOj spray 
even at the downstream end of the computation domain. Predicted temperatures are about 
3000°K, suggesting the presence of drops is unlikely. No mechanism has yet been suggested 
by these investigators or others for explaining L0 2 drops at these high gas temperatures. If the 
mixing and/or combustion is reduced in the model to allow for the presence of LO 2 , the 
predicted combustion efficiency would be quite low. For combustion experiments, the primary 
measurement of interest is temperature. However, temperature has not been reported, therefore, 
no modification to the CFD model was made. 

3.4 Tripropellant Shear Coaxial Injector Elements 

The third test case simulated with the homogeneous spray combustion model was a 
tripropellant (G0 2 /GH 2 /RP-1) shear coaxial injector flow, also examined by Santoro and his 
coworkers 17 . The schematic of the injector element is sketched in Figure 16. The hydrogen gas 
is assumed to be fully mixed with RP-1 and injected into the combustion chamber through the 
inner tube of 0.15 in. diameter. The mass flow rates of GH 2 and RP-1 are 0.015 lb/sec and 
0.15 lb/sec, respectively. The gaseous oxygen has a mass flow rate of 0.48 lb/ sec and is 
injected through the annulus, which has inner and outer diameters of 0.18 in. and 0.5 in., 
respectively. The combustion chamber is modeled for 7 in. downstream from the injector face. 
The chamber pressure is 500 psi and the temperature of all propellants is 540 °R. The 
combustion is represented with three global kinetics expressions, as shown in Table 3, to 
simulate RP-1 oxidation, soot formation, and soot oxidation process, where the RP-1 oxidation 
step produces H 2 and CO. The combustion was completed with a wet CO equilibrium model. 
Hence, ignition delays may arise only from the RP-1 oxidation step. 

The predicted temperature field and regions of soot formation are shown in Figs. 17 and 
18, respectively. The flow is predicted to be highly striated, but all of the RP-1 is vaporized 
in the first 3 in. of the combustor. Drops with a mean size of 20 fim were reported at a point 
6 in. from the injector face and 0.2 in. off of the centerline, even though the combustion 
efficiency is reported to be 97.5%. Drop size distributions and velocities were reported for this 
point; however, the total flowrate of RP-1 liquid cannot be estimated from these data. The 
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velocity predictions and measurements are in close agreement, even though insufficient test data 
are available to identify weak points in the analysis. The soot mass fraction predictions are very 
low. However, soot is predicted where it would be expected. The soot model can be easily 
tuned, if suitable test data were available. 



Figure 16. Schematic of Injector Element Used for RP-1/GH 2 /G0 2 tripropellants 
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3.5 Fastrack Injector Elements 

After the firing of early versions of the Fastrack motor, the liner near the injector face 
plate showed locally excessive ablation. The large ablation occurred in streaks in the axial 
direction with the number of streaks, as well as the distribution of streaks over the wall in the 
circumferential direction, being dependent on the design of the injector face plate 18 . 

The objective of the CFD analysis of the near-injector-flow was to determine a possible 
source of the observed ablation phenomenon. For this investigation, the flow within a region 
of the motor that was bounded by a section of the injector face plate and the wall of the 
combustion chamber was modeled in detail. The section of the injector face plate used in the 
CFD modeling involved a pair (doublet) of oxidizer holes, a pair (doublet) of fuel holes, and 
three holes for film-cooling flow. Recognizing the symmetric nature of the flow under 
consideration with respect to the azimuthal direction allowed the flow to be computed for only 
one half of the injector face plate element with symmetry boundary conditions imposed to 
account for the other half. The face plate section modeled in the CFD analysis can be identified 
by referring to Fig. 2. 

The orientation based on a coordinate system erected in the center of the injector face 
and some dimensions of the volume treated by the CFD analysis are shown in Figure 19. 



Figure 19: Volume Element Used for CFD 
Analysis. 
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The details of the face plate section bounding the volume plotted in Figure 19 are shown 
in Figure 20. As shown in Figure 20, the computational domain consists of one L0 2 injector 
hole, one RP-1 injector hole, and one and one half RP-1 coolant injector holes. 



Figure 20: Injector Holes in Modeled Face 
Plate Section. 


Although the cross-sections of the injector holes in the CFD modeling are essentially 
rectangular, the area of each hole is the same as the corresponding cross-sectional area of the 
(elliptical) hole in the physical case. 

The computational grid consisted of 31 planes in the axial direction, with each plane 
employing 41 nodes in the radial direction, and 21 nodes in the azimuthal direction. A non- 
uniform grid spacing was selected in the axial direction with the grid points being packed more 
closely towards the injector face, thus allowing for a better resolution of the flow in the 
immediate vicinity of the injector. The grid points in both, the radial and the azimuthal are 
clustered in the injector hole regions to accurately capture the incoming mass-flows. The three- 
dimensional computational grid is shown in Figure 21. 
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Figure 21: Computational Grid Involving 
31x41x21 Nodes. 


A summary of the initial conditions for the different mass-flows entering the control 
volume is given in Table 12. The axial component of the velocity is implied by the mass 
flowrate. The other two components of the total velocity vector were calculated based on the 
orientation of each hole. The angles used for the calculations were taken from drawings 
supplied by NASA. The resulting total velocity vectors for the injected mass-flows are plotted 
in Figure 22. 


Table 12. The Inlet Flow Conditions for the Fastrack Injector Face Simulation 



L0 2 Hole 

RP-1 Hole 

RP-1 Coolant 

Area per hole (in 2 ) 

0.007726 

0.002432 

0.001077 

Mass flow rate per hole (lb/s) 

0.4042 

0.1482 

0.0291 

Temperature (°R) 

190 

532 

532 

Density (lb/ft 3 ) 

67 

46.5 

46.5 

Pressure (psi) 

514 

514 

514 
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Figure 22: Total Velocity Vectors of Injected 
Mass-Flows. 


Two numerical simulations of the flowfield near Fastrack injector face were conducted 
to analyze the effect of injector cant angles on the chamber wall heating. The first case 
corresponds to that just described. In Table 13, the impingement and cant angles of the L0 2 and 
RP-1 streams for two simulated cases are shown. For both cases, a symmetry condition was 
employed for the two boundaries in the circumferential direction. The top boundary in the radial 
direction is a wall; while the bottom boundary in the radial direction was assumed to be 
symmetry plane. The chemistry model used for the tripropellant injector was also used for this 
case. The flow angle of case #1 is the original design. The numerical results of case #1, as 
shown in Figure 23, indicates that a strong cross-stream secondary flow structure and a 
streamwise flow recirculation was induced due to the radial velocity component resulting from 
the propellant cant angles as well as the slip boundary condition imposed on the bottom surface. 
The oxygen and RP-1 flows were well mixed downstream of the impingement surface, and then 
was carried as secondary flow back towards the impingement surface and toward the top wall. 
Also, the propellant mixture recirculated towards, and heated, the injector face, and followed 
the same secondary flow structure. As a result, the hot streak on the wall occurred between the 
RP-1 coolant holes and aligned with the propellant injector holes. Though the strength of the 
secondary flow could have been over-predicted due to the slip boundary condition assumption 
of the bottom surface, the induction of the secondary flow is reasonable because of the large 
opening between the impingement point and the top wall. 

To reduce the severity of the induced secondary flow, a second test case was simulated, 
for which the cant angle for both the oxidizer and fuel flow was eliminated. The elimination 
of the flow cant angle makes the assumed slip boundary condition of the bottom surface less 
critical. As anticipated, the numerical results, as shown in Figure 24 reveal that a much weaker 
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Table 13. The Impingement and Cant Angles of the Propellants for Both Cases 



Case #1 

Case #2 

L0 2 impingement angle 

30° 

30° 

L0 2 cant angle 

50° 

0° 

RP-1 impingement angle 

25° 

25° 

RP-1 cant angle 

27° 

0° 


secondary flow, compared to that in case #\ was induced. As a result, the wall temperature is 
much lower and more uniform. Meanwhile, the region of stoichiometric mixture is much 
thinner compared to that in case #\ and is formed between the L0 2 doublet and the RP-1 
doublet. This indicates that secondary flows enhance mixing, but a compromise must be made 
to reduce wall heating vs. optimizing performance. 

Further study to add the inner rows of injector holes and thus reduce the sensitivity to 
the assumption of the slip boundary condition is recommended. Also, a gas only solution for 
the boundaries of the detailed injector element should be considered. Other alternatives to this 
homogeneous spray analysis could be devised, however, this methodology is believed to yield 
the most realistic flowfield information for the computation required for the analysis. The 
computation time required for these analyses was long, but further optimization and some 
additional simplification could be realized. 
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4.0 CONCLUSIONS & RECOMMENDATIONS 

Conclusions drawn from this investigation are: 

(1) The FDNS-RFV code is a useful analysis tool for evaluating near injector-face 
flowfields of tripropellant and bipropellant liquid motor injector element designs. 

(2) Experiments performed to characterize injector element configurations should be 
conducted in closer concert with analytical studies to insure that they produce 
sufficient test data to serve as verification of the applicable flowfield analysis. 
In particular, local temperature measurements are essential to validate a 
combustion model. 

Recommendations for future use of this CFD design code are: 

(1) Additional cases be analyzed with the FDNS-RFV code to determine the best 
method of synthesizing this detailed analysis with appropriate analyses of 
bounding flows. The use of multiple rows of injector elements and of gas only 
(or gas/spray) coarse grid analyses to serve as boundary conditions for detailed 
injector element analyses should be investigated. 

(2) The extension of this homogeneous spray combustion model to account for, even 
in an approximate fashion, the presence of spray drops should be considered. For 
example, finite-rate vaporization within the homogeneous spray combustion code 
is a simple extension which could prove quite useful. 
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APPENDIX A 
RADIATION EFFECTS 


• Radiative Properties along a Line-of-Sight 

Since scattering is negligible, the equation of transfer for a line of sight, s, through an 
absorbing and emitting medium, which is in local thermodynamic equilibrium, is used to identify 
and define the radiative properties. 

(dl.lds) + *,/, = k,l b , (1) 


where I v is the spectral intensity, I bv is the Planck blackbody intensity, k, is the spectral 
absorption coefficient and v denotes wave number. Integration of equation (1) over all wave 
numbers and along the line of sight from the point s = L to the observer at s = 0 yields 

w L 

i '{°} = || W ex P 


k ds' 


ds dv 


(2) 


I {0} is the spectral intensity directed toward point 0 along line s. 


The development of this equation closely follows that presented by Buckius & Tien [1]. 
A medium containing more than one component is characterized by the spectral absorption 
coefficient that is the sum of the spectral absorption coefficients of each species. With the 
transmissivity defined as 


r = exp 


ft,*- 


(3) 


the transmissivity for a mixture of i components is 


} = n 

i 


(4) 


which is strictly valid only for monochromatic radiation although it is found to be true 
experimentally over finite spectral intervals [2]. The intensity in equation (2) becomes 


/ { 0 > = [ | /„{»} 


dT r {s) 

ds 


dsdv 


(5) 
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• Homogeneous Gray Model 

The absorption and emission of radiant energy by a mixture of gases and particles is a 
function of the concentration, temperature, and pathlength. The radiative transfer calculation 
for a line of sight is, therefore, dependent upon the variation of these parameters along the 
pathlength. This can result in very complicated computations for mixtures of many components. 
For simplification, the temperature and pressure are often taken as constants resulting in a 
homogeneous path approximation. Another major simplification in the analysis is achieved if 
it is assumed that one mean value of the absorption coefficient can account for all spectral 
variations. These two approximations yield the homogeneous gray model. 

For the absorption coefficient which is not a function of wave number or pathlength, 
equation (5) becomes 

I- {0} = (oT 4 hr) [ 1 - exp {- KL }] (6) 

where K is a total averaged absorption coefficient. The use of this homogeneous gray equation 
is dependent upon the definition of K. Experimental measurements are one method of 
determining the mean coefficient. For both soot and combustion gases K is a strong function 
of v , therefore this simplification is not acceptable. 

• Homogeneous Nongrav Model 

The intensity for a homogeneous path is obtained from equation (5) as 

X 

/-{0} = f l„.\\-T,{L)]dv (7) 


The spectral variation of the transmissivity for an infrared radiating gas is due to the vibrational- 
rotational bands. The spectral dependence can be evaluated for practical calculations by either 
narrow or wide band methods [3]. In particular, the exponential wide band model [4] has been 
successfully employed for various combustion gases. The narrow band models are much more 
computationally intensive [5], but more importantly parameters for these models have been 
established from data on atmospheric and subatmospheric pressure flames. The wide band 
models are based on higher pressure flame data. 

For a gaseous medium with a single vibrational-rotational band, the frequency 
dependence for a band is restricted to a small interval as compared to the blackbody spectrum. 
For this reason, 1^ is removed from the integral in equation (7) and taken constant at the 
wavelength of the band head or center depending upon the type of band to give 
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The use of equation (8) implies that the gas absorbs like a blackbody over the band frequency 
interval. This procedure yields the correct band absorbance, but cannot be conveniently 
averaged with the continuum soot absorbance in this frequency interval. Edwards [6] remedies 
this problem by defining a band transmissivity, r g , with the equation 

T g = d In A * Id In t h < 0.9 (9) 

where r H is the maximum optical depth at the band head. Hence, 

A « = (^gjX^uj-^Ajj* = wAjj* 

Further explanation of this process is shown in an example problem given by Edwards [4]. 

Use of wide band formulas for the integrated band intensity a s , the band width parameter 
a);,, the mean line width-to-spacing ratio j8 ij5 and the equivalent broadening pressure P ei in the 
transmissivity of equation (9) allows the total band absorptance to be determined for each gas 
band. When two or more gases have bands in the same spectral interval, an overlap correction 
must be introduced. With equation (8), the total band absorptances, and the overlap correction 
must be introduced. Data for using Edwards exponential wide band model taken from [3] are 
shown in Table 1. Variations of this model are given in [7], 

As mentioned, narrow band radiation data and models for combustion gases are also 
available. The SIRRM [8] data base is build on empirical data taken with modest spectral 
resolution. The experiments were conducted primarily at atmospheric and sub-atmospheric 
pressure to evaluate rocket plume radiation. Recently, studies have been made to model 
molecular radiation directly by modeling all of the lines in a band. This has resulted in the 
HTTRAN data base [9-12]. Phillips [9] made a comparison between the band predictions using 
the SIRRM and HITRAN data; his work is still in progress. Several investigators have used 
the individual line data to predict and measure radiation at high pressures and modestly high 
temperatures (about 800 °K) [13-16]. Such data are very interesting to compare to wide band 
predictions, but require too much computational time to be utilized directly for simulating 
combustor heating. Leckner [17] compares integrations of narrow band model data for H 2 0 and 
C0 2 to mean absorption coefficient measurements. 

To account for soot radiation, a medium containing only particles may be considered. 
Small particles present in a medium attenuate radiation continuously throughout the spectrum. 
For small particle absorption, of any size distribution, the governing equation for the absorption 
coefficient is given as [1, 8] 

k _ 3 Qab^ r _ 36t /, nk (io) 

" 2 2 r J ’ X [n 2 -k 2 +2] 2 +(2nk) 2 

where k, p is the particle absorption coefficient, r the particle radius, Q, bs the absorption 
efficiency, f, the volume fraction of particles, and m = n-ik the index of refraction. The actual 
wavelength variation indicated in equation (10) is dependent upon the variation of m for the 
specific particles. For soot particles that are contained in combustion systems, the optical 
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Correlation Parameters 

Oq 

[cm'^g/m 2 )] 

7o 

Wo 

[cm 1 ] 


Table 1. Wide band model correlation parameters for various gases 


Pressure Parameters 


Band Location 


(*») 


HjO m = 3, rj, = 3652 cm' 1 , tj 2 = 1595 cm' 1 , - 3756 cm ', g k = (1, 1, 1) 


71 nm ' i) c = 140 cm' 1 1 8.6 (T 0 /T) 05 + 0.5 44,205 

(Rotational) 

( 0 , 0 , 0 ) 


6.3 nm tj c = 1600 cm' 1 1 8.6 (T 0 /T) 0 5 + 0.5 41.2 

(0,1,0) 


2.7 jtm t] c = 3760 cm' 1 
( 0 , 2 , 0 ) 
(1,0,0) 
(0,0,1) 


0.14311 


0.09427 


1 8.6 OVT) 05 + 0.5 



1.87 ij c = 5350 cm' 1 1 8.6 (To/T) 0 5 + 1.5 

( 0 , 1 , 1 ) 


1.38 fim ij, = 7250 cm' 1 1 8.6 (Tq/T ) 05 + 1.5 

( 1 , 0 , 1 ) 


CO, m = 3, rj, = 1351 cm ', i\ 2 = 666 cm', 173 = 2396 cm ', g k = (1, 2, 1) 


15 fim ij c = 667 cm' 1 0.7 
(0,1,0) 


10.4 /im tj c = 960 cm' 1 0.8 
(- 1 , 0 , 1 ) 


9.4 /im rj c = 1060 cm' 1 0.8 
(0,-2, 1) 


4.3 /xm ij c = 2410 cm' 1 0.8 
( 0 , 0 , 1 ) 


2.7 71m tj c = 3660 cm' 1 0.65 
( 1 , 0 , 1 ) 


2.0 /tm i} c = 5200 cm' 1 0.65 
(2,0,1) 


CO m = 1, tj, = 2143 cm' 1 , g, = 1 


4.7 fim i} c = 2143 cm' 1 0.8 
( 1 ) 


2.35 pm i} c = 4260 cm' 1 0.8 
( 2 ) 





2.47xlO r ® 0.04017 


2.48x1(3* 0.11888 


0.24723 


0.13341 


0.39305 




0.07506 


0.16758 



0.13219 b,c 60.0 b 


0.08169 


0.11628 


0.06157 
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Footnotes for Table 1. Parameters for the EWB Model 

a. For the rotational band a = a 0 ex P ("9 (Tf/T) 0,5 ), 7 = 7 0 (T(/T) 0 5 

b. Combination of three bands, all but weak (0,2,0) band are fundamental bands, a 0 
= 25.9 cm 1 

c. Line overlap for overlapping bands from 


a 




P = Y p e = Y 0 ^ 


2® ip p 
T *0 e ’ * 


-H-(l + (b-1) ^ 

Po P 


T 0 = 100°K, po = 1 atm 


«(D = a 0 


Y*(T) 

V'Vo) 


|l - exp 

f m > 

-E u k 6 k 
k-l 7 

" Y(T) 


m > 


1 - exp 

-E^* 8 k 

, k-i ; 

" YCTo) 


T(T) 


fr (v k + g k + \ - i)! c -u t v t 
k-lvX (g k - DlvJ 

TT y' ^ V k + &k ~ ^‘ e ~°k v k 

k=iv7o (g k -D!v k ! 


$ CD = 


A* = 

To = 

A* 

A* 

A‘ 

A* 

A* 


m 00 

P E , 

k-1 v k -v 0k \ 

(v k + g k + 8 k - l)! e - Uk v k 

(g k - l)!v k ! 

m « 

ns 

k-1 v k =0 

(v k + g k + « k - D! c -u k v k 

(g k - D!v k ! 

= A/co 

= a X/co 

X = < 

= T 0 

for /3 


2 (r 0 /3) 05 - 0 
In {r 0 /3} +2-/3 


= T n 


In t c + 1 


for /3 < 1 1//3 < r 0 < 00 
for /3 < 1 and 0 < r 0 < 1 
for /3 < 1 and 1 < r c < 00 
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properties have been recently reviewed by Reid [18]. The review indicated that several 
dispersion models for m have been reported. Furthermore, although the dispersion models 
indicate significant differences in m the predicted volume IR absorption coefficients vary much 
less. These investigators chose the model of Lee and Tien [19] for coding to generate Figure 
1. These index of refraction values are recommended until better data become available. The 
data are believed to be reasonably accurate. 


Using this radiation model for soot, the rocket motor test data of Boynton [20] were re- 
analyzed. The soot density was determined by considering the soot to be the sole radiator at 
2.2 /tm. The re-analyzed data are shown in Figs. 2-4 for three nozzle configurations. The 
excellent fits of the test data show that the results are consistent with the Lee & Tien radiation 
model. The test report suggested that gaseous radiation was estimated and removed form the 
particle radiation, but details of this removal process were not specified. Additional analyses 
of these data are presented in the Phase II progress reports [21] that this Phase III study 
supports. 


Using equation (9) in the transmissivity given in equation (3) results in 


r, p {s] = exp 


3 Qabs 
2 2 r 


s 



( 11 ) 


For a homogeneous path, f, is a constant and the integral in equation (11) is f, s. This result 
gives the intensity for a particulate medium in terms of the particle’s properties, temperature, 
and pathlength. This can also be used to determine the temperature of the homogeneous path 
if f„ Q»b.> and L are known. 

If the medium is composed of both particles and gases, the overlapping spectra must be 
analyzed. Particulate emission and absorption is continuous and is essentially constant over the 
banded gas regions. Therefore, the particle transmissivity can be removed from integrals over 
gas bands and evaluated at the band center or head. The governing equation for the total 
intensity is 


x 

r i°} = ^h.v- T . P ) dv + 


( 12 ) 


where Ay accounts for the overlap of species with one another and is the transmissivity of 
the particles. The first term denotes the intensity resulting from particulates, and the second 
term denotes the gaseous contribution that is transmitted through the particles. This is the result 
for the intensity leaving a homogeneous path containing gases and particles. The applicability 
of this result is dependent upon the degree of homogeneity present in the system. 

Edwards [22] presented results of this type simulation for conditions typical of a rocket 
motor with pressure level as a parameter. The highest pressure which he considered was about 
950 psi. At this pressure the radiation is still not blackbody (for the soot level which was 
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Figure 2. 14 Foelsch Nozzle Radiation Measurements - High Expansion Ratio 
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Figure 4. Foelsch Nozzle Radiation Measurements - Low Expansion Ratio 
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assumed). This result is shown in Fig. 5, and it suggests that his radiation property estimation 
method might be useful for predicting radiation within rocket combustors. Had blackbody 
radiation been predicted at this and lower pressures, the method would have been suspect and 
further consideration been recommended. 

• Non-homogeneous non gray model 

The most complete method to determine the intensity from a general pathlength is to 
incorporate the variations in temperature and pressure. The analysis needed is, therefore, more 
complicated, and more detailed information is required. The resulting expressions are used for 
pathlengths that have large variations in temperature and pressure where the homogeneous 
nongray model is no longer valid. 

The general non-homogeneous intensity at s = 0 is given in equation (5). Similar to 
equation (12) for a homogeneous nongray path, the non-homogeneous nongray model for a 
mixture of particles and gases if given as 



L 

+ p ~ 
0 



0 


fds 


( 13 ) 




n°) = / p bv -j [1 - [s)}dsdv + /£ I bvij [^ pii A.]ds 

00 “ 0 V 


( 14 ) 


where the overlap correction is introduced into Ay. The terms of equation (14) are interpreted 
in a similar manner to the homogeneous result equation (12) except that all quantities are now 
pathlength dependent. This is the general governing expression for the non-homogeneous 
nongray model. 

An approximate method of averaging to replace non-homogeneous regions with scaled 
values of the wide band parameters for use in the homogeneous total band absorptance has been 
reported [23]. This work gives the average values as 

a - — f ads 1 


The use of any of the scaling techniques requires detailed pathlength information about the 
variation of the partial pressure and temperature. 

• Solution of the Uncoupled Radiative Transfer Equation 

The radiative transfer equation (RTE) is an integro-differential equation that expresses 
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Figure 5. Combustion Gas and Soot Radiation 
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the conservation of radiation. The concentration and temperature fields are needed to evaluate 
this equation. From a solution of the RTE a radiation heat flux vector (q,) can be obtained 
which enters directly into the energy equation. Usually q, is small enough that it can be 
neglected in the energy equation, such that the RTE can be solved in an uncoupled fashion after 
the remaining conservation equations are solved. This is most fortunate, because the RTE cannot 
be solved in general, even numerically. The closest approach is to use a Monte Carlo [24] 
simulation to obtain the RTE solution. It is not feasible to obtain a Monte Carlo solution at each 
grid point at each CFD time step to generate a truly coupled solution of the entire system of 
conservation equations. This dilemma has resulted in a multitude of approximate RTE solutions 
being used to represent the radiation in the expectation that some form of "coupled" solution can 
be obtained. It is absolutely essential that realistic line-of-sight (LOS) radiance calculations be 
made in order to determine the extent to which approximations can be accepted in RTE 
solutions. The radiation properties presented in the previous section are sufficiently accurate to 
made the LOS analyses using the RTE solution technology described in the following paragraph. 
The critical issue is the accuracy of the prediction of gaseous species and soot concentrations 
which characterize tripropellant combustor flowfields of interest. 

SECA has developed and/or collected and utilized several codes for use in predicting 
radiation. These codes are listed and described in Table 2. The Lee and Tien soot data are 
indicated as L&T soot in Table 2. The NOZZRAD [25], GASRAD [26], and SIRRM [8] codes 
are particularly suitable for evaluating LOS radiance characteristics. The NOZZRAD code can 
be used as a subroutine to provide LOS data for use in other radiation analyses such as IDA 
[27]. The IDA code is best method of offering a fully coupled radiation analysis, provided the 
LOS analyses justify the assumptions used in its development. The diffusion approximation code 
(DAC) [27] offers a very convenient method of predicting wall heating, again provided its use 
can be justified. 

The LOS analyses were planned as part of this study, but the reduced effort resulting 
from funding limitations and the fact that all of the other aspects of the study had to be 
performed first in order to obtain realistic concentration fields (especially of soot) prevented 
these analyses from being made. Upstream radiative heating of the flowfield near the injector 
face and of the injector face itself is extremely important for rocket motor design, but these 
flowfields are very difficult to predict due to geometrical complexity and the possible presence 
of significant amounts of liquid droplets. Hence, their prediction is a costly and uncertain 
endeavor. 

• Coupled Radiative Transfer Equation Solutions 

The extensive furnace modeling efforts of Lockwood and his colleagues [28] and flame 
modeling studies at the Naval Research Laboratory [29] and by Moss, et al [30] indicate that 
CFD/radiation analyses are feasible. However, none of these studies simulated conditions which 
are typical of rocket motor combustion, and all of the radiation analyses contain approximations 
which must be re-examined before they can be applied to rocket motor combustion chambers. 
The flame studies performed at the National Research Council of Canada [31] are very valuable 
for verifying concentration predictions of tripropellant combustion since hydrogen enriched 
hydrocarbon flames were investigated by radiation experiments using a CARS system. 


A-13 



SEC A-FR-97 -04 


Table 2. RTE Solvers 

Code 

Methodology 

Geometry for 
Non- 

Homogeneous 

Fields 

Optical 

Properties 

■ 

NOZZRAD 

emission /absorption/ 
scattering 

i. LOS 

ii. Slabs 

L&T 

soot/Al 2 0 3 
particles wide 
band model 

25 

SIRRM 

emission/absorption/ 

scattering 

i. LOS 

ii. 2-& 6 Flux 

old soot data, 
A1 2 0 3 , narrow 
band model 

8 

GASRAD 

emission/absorption 

LOS from 
cylindrical field 

old soot data, 
narrow band 
model 

26 

IDA 

emission/absorption/ 
scattering w/ordinary 
differential 
approximation 

Cylindrical 

coordinates 

A1 2 0 3 , wide 
band model 

27 

DAC 

emission/absorption/ 
scattering diffusion 
limit 

Cylindrical 

coordinates 

ai 2 o 3 L&T 
soot, wide 
band model 

27 

REMCAR 

emission/absorption/ 
scattering, Monte 
Carlo 

3-dimensional 

A1 2 0 3 , old 
soot data, 
narrow band 
model 

24 
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Nomenclature for Appendix A 

Ay total band absorptance for jth band of the ith species 
c concentration of particles 

C 2 Planck’s second constant 

f v volume fraction of soot particles 

I intensity 

I" intensity in the negative s direction 

k complex part of refractive index 

k, absorption coefficient 

K average absorption coefficient 

L pathlength 

m complex index of refraction equal to n-ik 

n real part of refractive index 

P pressure 

P ei equivalent broadening pressure of the ith species 
Q, bs particle absorption efficiency 

r particle radius 

s distance along line of sight 

T temperature 

X mass pathlength 

Greek symbols 

a;, integrated band intensity 

x times the mean line width-to-spacing ratio 
e emissivity 

X wavelength 

v wave number 

p solid phase density of soot 

Pi density of ith absorber 

r transmissivity 

«jj band width parameter 

Subscripts 
b blackbody 

i species 

j band of specific species 

1 lower wave number limit 

p particle 

sc scaled parameter 

u upper wave number limit 

v wave number 

Superscripts 

(overbars) average values 
{ } denote functionality 
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